#
# Copyright (c) 2014 Luciano Xumerle. All rights reserved.
#
# This file is part of gephinfo.
#
# gephinfo is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
#

package bedInterval;

require 5.003;
use warnings;
use strict;
use utf8;

=head2 new( $bedFile )

Read intervals from bed e return the hash:

 $BED{CHROM}->{START}->{STOP}=1

and the hash:

 $BED{CHROM}->{STOP}->{START}=1

=cut

sub new
{
    my $class = shift;
    my $file  = shift;

    my $self = {};

    if ( defined $file && -f $file )
    {
        print STDERR qq|Start to import BED...\n|;
        my $intervals = 0;
        open( BEDFILE, $file ) || die;
        while (<BEDFILE>)
        {
            chomp;
            my @a = split /\t/;
            if ( defined $a[2] && $a[2] =~ m/^\d+$/ )
            {
                $self->{'STARTS'}->{ $a[0] }->{ ( $a[1] + 1 ) }->{ $a[2] } = 1;
                $self->{'STOPS'}->{ $a[0] }->{ $a[2] }->{ ( $a[1] + 1 ) } = 1;
                $intervals++;
            }
        }
        print STDERR qq|$intervals imported from BED file!\n|;
        close BEDFILE;
    }

    bless( $self, $class );
    return ($self);
}

=head2 isPresent( $chrom, $start, $stop )

Returns true (1) if the interval exists. Else false (0).

=cut

sub isPresent
{
    my $self  = shift;
    my $chrom = shift;
    my $start = shift;
    my $stop  = shift;

    return 1
      if ( defined $self->{'STARTS'}->{$chrom}->{$start}->{$stop}
        && $self->{'STARTS'}->{$chrom}->{$start}->{$stop} == 1 );
    return 0;
}

=head2 getStop( $chrom, $start )

Returns true smallest end base for the given start position.
Return -1 if the start base doesn't exist.

=cut

sub getStop
{
    my $self  = shift;
    my $chrom = shift;
    my $start = shift;

    my @a = keys %{ $self->{'STARTS'}->{$chrom}->{$start} };
    if ( defined $a[0] )
    {
        @a = sort { $a <=> $b } @a if ( $#a > 0 );
        return $a[0];
    }
    return -1;
}

=head2 getStart( $chrom, $stop )

Returns true greater start base for the given stop position.
Return -1 if the stop base doesn't exist.

=cut

sub getStart
{
    my $self  = shift;
    my $chrom = shift;
    my $stop  = shift;

    my @a = keys %{ $self->{'STOPS'}->{$chrom}->{$stop} };
    if ( defined $a[0] )
    {
        @a = sort { $b <=> $a } @a if ( $#a > 0 );
        return $a[0];
    }
    return -1;
}

=head2 getIntervalsMinorThan( $length )

Returns the list of bed intervals with length less than $length.
Each bed inteval is reportes as list like: [ chrom, start, stop ].

=cut

sub getIntervalsMinorThan
{
    my $self = shift;
    my $size = shift;
    my $res  = [];

    foreach my $ch ( keys %{ $self->{'STARTS'} } )
    {
        foreach my $sta ( keys %{ $self->{'STARTS'}->{$ch} } )
        {
            foreach my $end ( keys %{ $self->{'STARTS'}->{$ch}->{$sta} } )
            {
                push @$res, [ $ch, $sta, $end ]
                  if ( $end - $sta + 1 < $size );
            }
        }
    }
    return $res;
}

=head2 getSequenceFromReference( $chrom_bin_dir, $chr, $start, $stop )

A static method to retrieve a subsequence from the reference sequence.
The method uses the chromosome in bin format: each file is named as
chrN.fa.bin and contains the sequence chromosome on a sigle line
without spaces charachters.

=cut

sub getSequenceFromReference
{
    my $chrom_dir = shift;
    my $chr       = shift;
    my $pos       = shift;
    my $stop      = shift;
    $stop = $pos if ( !defined $stop );

    if ( $chr !~ m/^chr/ )
    {
        $chr = 'chr' . $chr;
        $chr =~ s/T//;
        $chr =~ s/chrXY/chrX/;
    }

    $chrom_dir =~ s|/$||;
    my $res = '';

    open( FILEFILE, $chrom_dir . '/' . $chr . '.fa.bin' ) || return '';

    for ( my $i = $pos ; $i <= $stop ; $i++ )
    {
        seek FILEFILE, $i - 1, 0;
        $res .= uc( getc FILEFILE );
    }
    close FILEFILE;
    return $res;
}

=back

=head1 AUTHOR

Luciano Xumerle (luciano.xumerle@gmail.com)
